Energy radiation of moving cracks 
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The energy radiated by moving cracks in a discrete background is analyzed. The energy flow 
through a given surface is expressed in terms of a generalized Poynting vector. The velocity of the 
crack is determined by the radiation by the crack tip. The radiation becomes more isotropic as the 
crack velocity approaches the instability threshold. 



I. INTRODUCTION 

The dynamics of cracks in brittle materials are being 
extensively studied |l],|2| , and a wealth of instabilities and 
patterns have been observed as a function of control pa- 
rameters such as the applied strain , or thermal gra- 
dients H . The theoretical analysis of moving cracks was 
initiated long ago J§-[ll[], with the study of exact solu- 
tions for cracks moving at constant velocity. These stud- 
ies have been extended to a variety of different situations 
[Ejyni. Alternatively, analytical approximations to the 
leading instabilities of a moving tip have been proposed 
& 

The simplest discrete model which captures the main 
features of cracks in brittle materials is a a lattice with 
central forces (springs) between nearest neighbors, whose 
bonds lose the restoring force above a given threshold 
fl5f (for extensions see also [^6|). This model, or sim- 
plifications of it which leave out the vectorial nature of 
the strain field, has been extensively used in modelling 
moving cracks HJl^jOHDl > although models which deal 
with the microscopic structure of the system are also be- 
ing considered ]2(i|-p3| . Alternatively, various continuum 
models, which describe the fractured zone in terms of 
additional fields, have been proposed p^ , p5[ . 

Discrete and continuum models of cracks differ in a 
variety of features. It is known that the discrete models 
used so far cannot describe a fracture zone at scales other 
than the size of the lattice cell in the calculations Jl8],^6| , 
although, even for a canonical material such as PMMA, 
the fracture zone has a dimension much larger than the 
size of its molecular building blocks J27|] . 

Another important difference between discrete and 
continuum models is the existence of radiation from the 
tip of the moving crack, due to the existence of periodic 
modulations in the velocity in the presence of an underly- 
ing lattice. In this sense, a lattice model for cracks is the 
simplest example where radiation due to the scattering 
of elastic waves by deviations from perfect homogeneity 
can be studied. These processes have been observed in 
experiments [|s|j28|~|3C|] , and it has been argued that they 
are responsabible for some of the crack instabilities BG] . 



In the present work, we study the energy radiated by 
a crack moving at constant velocity in a discrete lattice. 
We use a generalization of the scheme discussed in [0. 
The general method used is explained next. Section III 
presents the main features of the results. The physical 
implications of the results is discussed in section IV. 

The problem of sound emission by moving cracks has 
been addressed, within a different scheme, by Inso- 
far as the two approaches can be compared, the results 
are compatible. Finally, radiation of moving cracks along 
the edge of the crack can be important in understanding 
the roughness of the crack surface |52| . We will focus on 
the radiation along the crack surface, and into the bulk of 
the sample. Experiments |m| and simulations suggest 
that this type of radiation can play a role in the observed 
instabilities of the crack tip. 



II. THE METHOD 

We study discrete models of elastic lattices in two di- 
mensional stripes, as discussed in |l^,^8|. The underlying 
lattice is hexagonal, with nearest neighbor forces. Bonds 
break when their elongation exceed a given threshold, 
Uth, and under a constant strain at the edges, which, 
scaled to the width of the stripe, we write as u . We 
study models with and without dissipation in the dy- 
namics of the nodes. Results depend on the ratio Uth/uo 
(see section III for further details). 



A. Energy considerations 

In the absence of dissipation, the total kinetic plus elas- 
tic energy must be conserved. In a continuum model, in 
the absence of radiation, energy conservation leads to a 
global and to a local constraint, for cracks moving at 
constant speed, v. 

i) In the absence of radiation, the region well behind 
the crack tip has relaxed to equilibrium, while the re- 
gion ahead of it is under the applied strain. The relaxed 
region grows at the expense of the region under strain, 
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at constant rate oc u^vW where W is the width of the 
stripe Jl8| . Energy is transferred to the crack at this rate. 
As the energy stored in the crack grows at rate oc u 2 h v, 
the crack can only propagate (without radiation) for a 
fixed value of uth/uo. Note that continuum solutions for 
radiationless cracks moving at constant speed Q do not 
specify a parameter equivalent to u t h, so that they do 
not conflict with energy conservation. 

ii) The only position at which elastic energy is used to 
increase the size of the crack is the crack tip. Thus, at 
the crack tip the flux of elastic energy should be equal to 
the energy invested in enlarging the crack jl],[33| . In con- 
tinuum models, this local constraint leads to an equation 
of the type: 



T = A(v)G 



(1) 



where T is the crack energy per unit length (T oc uth 
in our lattice model), G is proportional to the stress in- 
tensity factor at the crack tip (G oc Uq), and A(v) is a 
universal function which goes from 1 at v — to at 
v = vr where vr is the Rayleigh speed. 

In lattice models, the energy arguments have to be 
modified because of the presence of radiation of elastic 
waves. If we assume that the difference between contin- 
uum and lattice models is small, we can use the pertur- 
bative scheme discussed in j34|-|38ll. The crack tip veloc- 
ity undergoes oscillations at frequency to = v/a, where 
a is the lattice constant, and amplitude /. In order to 
estimate the energy radiated from the tip, we have to 
extend the perturbative expansion to second order. We 
will not attempt here to calculate this expansion rigor- 
ously. However, from the knowledge of the leading term 
p6|p8[ , we can infer that the radiation due to a pertur- 
bation of frequency uj should go as Buj 2 , where B is a 
positive constant. The global constraint i) implies that 
f 2 (v/l) 2 oc 1 — k(uth/uo) 2 , where I is the crack length 
and A: is a constant. The radiation also has to be taken 
into account when balancing the energy absorbed at the 
crack tip, ii) above. The corrections to eq. ([!]) lead to a 
condition of the type: 
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where c and d are constants, and A is proportional to 
the inverse of the function A defined in eq.(Q), the magni- 
tude in brackets being its argument. Thus, the existence 
of radiation in discrete models allows for the existence of 
a continuum of solutions v(u /u t h) for a given range of 
u 1 P^ , p^ , p^7| -p9[ . Note that the existence of solutions 
which do not violate energy conservation does not imply 
that these solutions are stable. Full dynamical simula- 
tions of lattice models |17 19 suggest that inertial cracks 
(without dissipation) accelerate until they reach speeds 



comparable to those predicted by the Yoffe criterion 
and then bifurcate. 



B. Energy flux: continuum elasticity 

In the following, we will reformulate the concepts dis- 
cussed in in order to make them more amenable for 
extensions to lattice models, dicussed in the next subsec- 
tion. 

We describe an elastic medium in terms of the energy 
& 



H = Hh 



Helas 



H kin - jd r-l-^- 



(3) 



where D is the spatial dimension, p is the mass density, 
A and u are Lame coefficients, u(r) denotes the displace- 
ments at position r, and the « tJ 's define the strain tensor: 



if dui duj 
2 \drj dn 



(4) 



The equations of motion satisfied by u(r) can be written 
as: 



d 



dt 2 



(•5) 



where Uji = dTi. e i as tic/9uij is the stress tensor. 

The time derivative of the total energy Eq within a 
region Q is: 



dE r , 
dt 
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(6) 



so that the vector P(r) with components Pj — 
X]j (Jijdui/dt plays the same role as the Poynting vector 
in electrodynamics. The energy flux through an element 
of area dS is given by PdS. Note, however, that, un- 
like in clcctromagnetism, the equations of elasticity have 
not Lorentz invariance (there are two sound velocities), 
and it is not possible to define a four vector combining 
P and the energy density. The energy transferred to the 
outside of this region remains defined as the flux of the 
vector P through the surface bounding Q. In the pres- 
ence of dissipation, we still use P as defined in Eq. (|^) 
in the understanding that what viscosity does is to trig- 
ger the partial absorbtion of the radiated enegy whithout 
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changing the direction in which it is emitted. The vector 
P will be our starting point in the study of the energy 
flux of a moving crack. 

C. Energy flux: lattice model 

We will compute numerically the radiation of energy 
in a discrete model, defined as a hexagonal two dimen- 
sional lattice with nearest neighbor forces |^5| , p7|Jl8| . The 
energy is given by the sum of a kinetic term, associated 
to the velocities of the nodes, and an elastic term, due 
to the deformation of the bonds. The variation of the 
elastic energy of a given bond with time can be written 
as: 

^-klv-WW^MM (7 ) 

where k is the force constant, and is a unit vector 
in the direction of the bond. We distribute this energy 
among the two nodes connected by the bond, so that we 
can write the total elastic energy within a given region 
as a sum of the contributions of the nodes within that 
region. The variation in the kinetic energy at node i is: 

3 

The variation of the total energy within a given region 
is calculated by summing over all bonds within that re- 
gion. The kinetic and elastic contributions for all bonds 
outside the edge of the region cancel. We are left with 
surface terms only, as in the continuum model described 
earlier. The surface contributions can be written as a 
sum of terms associated to the bonds which connect a 
node within the region under study and a node outside. 
Thus, a surface which includes a given node and has a 
given orientation leads to an energy flux across it which 
can be calculated from a weighted sum of the positions 
and velocities of the bonds which connect that node to 
its neighbors. As we can associate to each surface orien- 
tation an energy flux, we can define the lattice Poynting 
vector, in analogy to the analysis done for the continuum 
model. We will use this discrete Poynting vector in the 
discussion of the energy dissipation of a moving crack 
below. 

III. RESULTS 

The discrete equations of motion in a two dimensional 
lattice of a given size are integrated numerically as dis- 
cussed in detail elsewhere [M . The lattice is maintained 
under constant load at the edges. In order to obtain 
cracks moving at constant velocities, a notch is induced 
at one side, which is gradually enlarged, along a straight 



line, until the stress buildup leads to the spontaneous 
propagation of the crack. The crack position, as function 
of time, is shown in Figure [I], for two different applied 
strains. The calculations show that the crack propagates 
freely at a constant velocity in the steady state. Our 
method for the calculation of the properties of cracks 
moving at constant speeds should lead to the same re- 
sults as other techniques. 




200 400 

time 



FIG. 1. Crack tip displacement versus time for cracks un- 
der two different applied strains (upper thick line: uo = 0.08; 
lower thick line: uo = 0.02) and zero viscosity. The arrow 
indicates the position of the notch beyond which the stresses 
at the crack tip exceed the threshold stress, and the crack 
propagates freely. The thin line is for uo = 0.08, but allowing 
for the branching instability. 

Instabilities are avoided by allowing only the bonds 
directly ahead of the crack to break. In other words, 
we force the crack to propagate straightly (with no 
branching). The simulations are performed in systems 
of 400 x 120 lattice sites, where we have checked that fi- 
nite size effects on the steady state velocity are less than 
1 percent. 

A. Crack velocity 

Figure || shows the steady state velocity v as a func- 
tion of the applied strain uq, for two different values of 
the viscosity (77 = and 77 = 0.8 in our units). The crack 
velocity increases monotonically with uq and asymptoti- 
cally tends to its limiting value = 0.571, the Rayleigh 
velocity in units where the force constant k = 1 and the 
mass per site m = 1 [fill . Due to lattice trapping, there 
is a minimum allowed uq whose value is roughly inde- 
pendent of 77 HQ, which in turn leads to a minimum 
crack speed which depends strongly on rj. The arrow 
marks the instability that would occur if the crack were 
not constrained to move on a straight line. 
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FIG. 2. Crack velocity versus external strain, for r\ = (up- 
per curve) and r\ = 0.8 (lower curve). The arrow on the right 
indicates the Rayleigh velocity vr. The vertical arrows mark 
the (avoided) branching instability (see text). The threshold 
for breaking is u±h = 0.1 



B. Elastic energy and hoop stress 




FIG. 3. Density of elastic energy (upper panel) and hoop 
stress (lower panel) for an inertial crack (77 = 0) moving under 
an applied strain Uo = 0.02, below the branching instability. 
The crack moves from left to right, and the tip is located at 
a point of coordinates 170, 60. 



Figure || shows snapshots of the density of elastic en- 
ergy and the hoop stress at a given time to-, for a steady 
inertial crack (r) = 0) moving at a velocity below the 
branching threshold (the applied strain is Uq = 0.02, cf. 
figure g). 

The density of elastic energy (FigM upper panel) has 
a sharp peak at the crack tip. In the near region (a few 
lattice spacings away from the tip), we see that the distri- 
bution of elastic energy is very anisotropic: it is sizeable 
in the direction perpendicular to the crack motion, where 
it decays smoothly with the distance, and all along the 
crack, where it has an oscillating behavior. This behav- 
ior is reminiscent of the Rayleigh waves which propagate 
on the crack surface (see section III C| below). At larger 
distances (of the order of the linear dimensions of the 
system), the elastic energy is smoother and has a broad 
maximum ahead of the tip, around a given angle of the 
order of 9 « 7r/3 from the crack direction. We cannot be 
conclusive about this maximum being intrinsic in nature, 
or rather being related to the symmetry of the underlying 
triangular lattice (see []l9| for a more detailed discussion 
of this point). 

The hoop stress (Fig. |^, lower panel) shows a very 
similar behavior, with strong oscillations all along the 
crack, and maxima perpendicular to the crack motion, 
the maximum shifting from 9 rs 7r/2 to 9 « 7r/3 with 
increasing distance from the tip. 
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FIG. 4. Density of elastic energy (upper panel) and hoop 
stress (lower panel) for an inertial crack (77 = 0) moving un- 
der an applied strain uq = 0.08, well above the branching 
instability. 

Figure ^ is the same as Fig. ||, but for a crack mov- 
ing at a velocity well above the branching threshold (the 
applied strain is u = 0.08). We notice that the distribu- 
tion of elastic energy and hoop stress has changed quali- 
tatively: the bulk features in the direction perpendicular 
to the crack motion now dominate over the oscillating 
part along the crack. The latter decay more rapidly and 
eventually disappear far behind the tip. 




FIG. 5. Density of elastic energy (upper panel) and hoop 
stress (lower panel) for a dissipative crack (77 = 0.8) moving 
under an applied strain uq = 0.08. 

The elastic energy and hoop stress corresponding to a 
dissipative crack (77 = 0.8, uo = 0.08) are shown in fig. 
pL Although the overall characteristics are similar to the 
inertial case, with maxima at the tip and in the direc- 
tion transverse to the crack, the distribution of stresses 
is much smoother. Moreover, the oscillations associated 
with Rayleigh waves along the crack are washed out by 
viscosity, being replaced by a single broad maximum be- 
hind the tip. 



C. Radiation 
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FIG. 6. Poynting vector field representing the radiation 
propagating in the vicinity of the crack tip. Upper panel: 
slow inertial crack (same parameters as fig. |^); center panel: 
fast inertial crack (same as fig. ^); lower panel: dissipative 
crack (same as fig. rjl). 

The above results can be better understood by analyz- 
ing the Poynting vector field, which represents the flux of 
energy being radiated at a given point in the system. As 
was stated in the introduction, emission of sound waves is 
expected since the crack tip moves in a discrete medium, 
therefore acting as a source of radiation at a frequency 
lu = v/a, the ratio of the crack speed to the lattice spac- 
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ing. Moreover, one expects a net flux of energy in the 
direction opposite to the crack motion, corresponding to 
the elastic energy released from the region ahead of the 
tip, which allows the crack to move. 

As can be seen in the first panel of figure ||, at such 
moderate crack speeds most of the energy is radiated in 
the form of Rayleigh waves propagating backwards along 
the crack, with a wavelength comparable with (but not 
equal to) the lattice spacing a. Despite the fact that 
1] = 0, such waves are seen to decay at long distances 
behind the tip (they decay into bulk waves, the oscil- 
lating bonds on the crack surface acting themselves as 
sources of radiation). In addition, there is also a weaker 
emission of bulk waves from the tip, responsible for the 
observed maximum in the direction perpendicular to the 
crack motion. 

At high crack speeds, on the other hand (cf. center 
panel in fig. ^J), it is the bulk radiation which dominates 
the emission pattern. Moreover, shadow images of the 
near-field appear behind the tip (the strongest one being 
at around x = 161). 

In the case of viscous cracks (lower panel in figure ||) , 
the emission pattern is entirely dominated by bulk waves, 
and Rayleigh oscillations disappear in agreement with the 
reults of fig. |[ 
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FIG. 7. The component of the Poynting vector along an 
inertial crack, as a function of the coordinate x, normalized 
to the total energy flowing through the sample, oc UqV W (see 
text). From top to bottom: uo = 0.02, 0.04, 0.08. The vertical 
arrow marks the tip position. The curves are shifted by a 
vertical offset for clarity purposes. 

In order to analyze the behavior of Rayleigh waves 
along inertial cracks, we plot in figure ^ the component 
of the Poynting vector P x parallel to the crack direction, 
at the surface of the crack, for different values of the ap- 
plied strain. The data are normalized to the difference in 
mechanical energy of a line far ahead from the crack, and 
a line far behind (this energy, which scales as UqvW, is 
transferred to the crack in the fracture process) . The fig- 
ure clearly shows that the wavelength of surface waves as 



well as their decay rate increase with the crack velocity. 

IV. CONCLUSIONS 

We have analyzed the nature and influence of radiation 
in the propagation of cracks in discrete systems. For the 
lattice and force models that we have studied, we find: 

i) Cracks in lattice models radiate energy, even when 
the average velocity is constant and they move along a 
straight line. This can be understood by assuming that 
the crack tip undergoes oscillations at frequencies nv/a, 
where v is the velocity of the crack and a is the lattice 
constant. 

ii) Radiation allows for the existence of a continuum 
of solutions of moving cracks at constant velocity. The 
balance of static elastic and crack energy is compensated 
by the radiation from the crack tip. 

hi) At low velocities, most of the radiation is in 
Rayleigh waves along the surface of the crack. At ve- 
locities comparable to the Rayleigh velocity, a significant 
fraction of the radiated energy is in bulk waves with a 
more isotropic distribution. 

iv) Viscosity allows for a faster exchange of the elastic 
energy stored ahead of the crack tip into other forms of 
energy. This can help to explain the increased stability 
of straight cracks in the presence of viscosity. 

Among the questions which remain unsolved is the re- 
lation of the radiation to the instabilities of the crack tip. 
Our results suggest that inertial cracks accelerate along 
a straight line, until they attain speeds compatible with 
Yoffe's criterion |1. On the other hand, the radiation of 
the crack tip becomes more isotropic at high velocities. 
It is unclear whether the continuum approach suffices to 
understand the instability observed in dynamical simula- 
tions of discrete models, or if the radiation from the tip 
of the crack plays a role in the instability. Note that the 
calculated instability occurs at higher velocities than the 
instabilities observed experimentally. 
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